Avalanche Size Scaling in Sheared Three-Dimensional Amorphous Solid 
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We have studied the statistics of plastic rearrangement events in a simulated amorphous solid at T = 0. 
Events are characterized by the energy release and the "slip volume", the product of plastic strain and system 
volume. Their distributions for a given system size L appear to be exponential, but a characteristic event size 
cannot be inferred, because the mean values of these quantities increase as L°' with a ~ 3/2. In contrast to 
results obtained in 2D models, we do not see simply connected avalanches. The exponent suggests a fractal 
shape of the avalanches, which is also evidenced by the mean fractal dimension and participation ratio. 
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Athermal, or low temperature, plastic deformation of amor- 
phous solids exhibits intermittent stress fluctuations and shear 
localization, in materials as diverse as metallic glasses 111], 
granular materials [2], foams f3\ and glassy polymers De- 
tailed knowledge of plastic deformation mechanisms in amor- 
phous solids, and their connection to macroscopic flow prop- 
erties, however, remains elusive: While in crystals the disloca- 
tion provides a well defined starting point for estimates of flow 
stress, a-f, in glasses there is no such easily characterizable 
defect. The traditional picture of deformation in amorphous 
solids-pioneered by Argon jstl and coworkers-is that plastic- 
ity involves collections of 'relaxation centers' [6] or 'shear- 
transformation zones'(STZs) |7] which operate as localized 
centers of deformation. This picture is supported by simula- 
tions of deformation in amorphous met als l|7i Isl Isl [lOl lnl ll2[] 
and observations of localized events Il3|] 7 

Mean-field theories of plasticity iHItIHi rely on this view- 
point, with the additional assumption that STZs operate some- 
what independently. But the detailed nature of correlations 
between shear transformations is a subtle issue. Elementary 
shear transformations should give rise to long-range elastic 
displacement fields, in analogy with the transformation of el- 
liptic Eshelby inclusions, in particular, al/r^ stress field due 
to a compact source. Models incorporating such interactions 
ifTsi [l6] exhibit localization of deformation in patterns remi- 
niscent of shear bands ||l7|] . suggesting that long-range elastic 
interactions may play an important role in the plastic response. 
The occurrence of shear bands in metallic glasses has been a 
major obstacle in the' development of these materials for en- 
gineering applications lllSfl . 

Even in carefully prepared samples (free of fracture- 
producing flaws), experimental observation of plastic defor- 
mation is often hindered by localization. In numerical simula- 
tions, however, it is possible to preserve translation-invariance 
and access statistical properties of plasticity in steady state. 
Athermal, quasi-static deformation allows further simplifica- 
tion of the underlying dynamics: Lacks has shown that in po- 
tential energy landscape (PEL), plastic deformation involves 
the destabilization of local minima along a single zero mode 
The PEL point of view brought hope that elemen- 



tary shear transformations could be identified with elementary 
transitions between minima in the PEL, as often held by STZ 
theories. This notion has been challenged, however, in recent 
simulations in two and three dimensions (2D and 3D). These 
show that individual plastic events present multiple substruc- 
tures, which are more compact and localized [|2lL i22ll . Mal- 
oney and Lemaitre [21] found that when visualized accord- 
ing to active atoms, plastic events-transitions between local 
minima-tended to be localized in one dimension but spread 
out along the other. This behavior led to an apparent scaling 
where the energy released in an event scaled as L, the linear 
system size. 



It is essential to know whether the 2D results of Ref. 02111 
should transfer to 3D. This is not obvious because the cor- 
relations that lead to an event taking place over an extended 
region may depend crucially on the power-law dependence of 
the elastic Green's function which is weaker in 3D. In this 
Letter, we report 3D simulations of a realistic model metal- 
lic glass, undergoing athermal quasi-static shear deformation. 
The main novelties of our work are: (i) the dimensionality; 
(ii) our use of realistic interactions potentials; (iii) our use of 
fractal analysis to characterize the geometry of avalanches in 
3D. Our main results are (1) we observe a scaling of event 
sizes with exponent close to 3/2; (2) visualization of typical 
large events indicates that the avalanches are no longer "sim- 
ply connected", partially localized avalanches, but rather are 
spread throughout the simulation box, with an apparent fractal 
shape, with mean fractal dimension close to 3/2. We note that 
this scaling behavior we differs from that observed in stud- 
ies of crackling noise in magnets, dislocation avalanches in 
single-crystal plasticity [23] or other systems that are charac- 
terized by critical behavior which leads to power law distribu- 
tions of event sizes Il24ll . For such systems finite-size effects 
only influence the large-avalanche tail of the power law dis- 
tribution. In our case the delocalized nature of the avalanches 
leads to a situation where the whole distribution scales with 
system size. 

The simulated material is Mgo.85Cuo.15, which is the ogti- 
mal glass-forming composition for the Mg-Cu system ll25ll . 
This system is interesting because the addition of a small 



amount of Y makes it a bulk metallic glass with high strength 
and low weight [261 . The interatomic potential is the effec- 
tive medium theory 12711 . fitted to properties of the pure ele- 
ments and intermetallic compounds obtained from experiment 
and density functional theory calculations. The configurations 
were created by cooling from a liquid state above the melt- 
ing temperature down to T=0, using constant temperature and 
pressure molecular dynamics. Details of the potential and of 
the cooling process may be found in 112 811 : the cooling rate 
for the systems studied here was about lO^^Ks^^. Periodic 
boundary conditions were used both in cooling and in the 
deformation simulations described below. Five system sizes 
were studied, containing 864, 2048, 4000, 8788 and 16384 
atoms, with L = 26A, 35A, 44A, 57A and 70A, respectively. 
Ten independent configurations were produced for all sizes, 
except only one 16384-atom system. 

The systems were deformed in pure shear by a varia- 
tion of the standard procedure of straining the entire sys- 
tem homogeneously in small increments followed by energy- 
minimization. In the so-called quasi-static limit, the system 
continuously follows deformation-induced changes in local 
minima 1 10, 19]; this protocol is meant to capture the asymp- 
totic trajectory in this limit. During homogeneous strains, as 
well as relaxing the atomic positions, components of strain 
apart from the one being controlled were also relaxed. In par- 
ticular this means that the hydrostatic pressure was always 
zero. To observe scaling it is important to be close to the 
quasi-static limit and therefore to have a strict tolerance for 
minimization 1,2 LI : we used lO^^eV/A for the maximum force 
and 10-^eV/A3 = 0.16 MPa for the maximum (relaxable) 
stress. The homogeneous strain was applied by multiplying 
the simulation box vectors by a shear strain matrix at each 
step, with a strain size of Ae — 0.0005. Here e is an off- 
diagonal component of the strain, not the engineering strain. 
The total strain at any point in the deformation history is the 
number of steps times 0.0005. The total amount of deforma- 
tion imposed was ^^100%. The use of periodic boundaries and 
the small sizes prevent any kind of macroscopic localization 
being observed even at such large strains. 

Examples of stress-strain curves are shown in Fig. [T] The 
behavior is very similar to that observed in other simulations 
m 2D and 3D ITS H 12111 : an initial linear elastic regime, 
followed by the onset of plasticity manifested as abrupt drops. 
After about 30% strain a steady state has been reached (note 
in particular the energy of the 16384-atom system). The stress 
averaged over the steady-state part of the curve, ct/, is about 
280 MPa, independent of L. 

The stress- and energy-drops define distinct plastic 
"events". These do not necessarily correspond to the elemen- 
tary units of plastic deformation, as mentioned above; a com- 
plete event is an avalanche of sub-events. Defining subevents 
in a useful way is problematic in practice, hence in order to 
stick with meaningful and well-defined quantities, we study 
complete events. 

We have analyzed the events in the stress and energy curves 
by assigning to each event two quantities. The energy and 
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FIG. 1: Stress and energy versus strain for 2048- (left) and 16384- 
atom (right) systems. 



stress drops are defined relative to the values they would be 
expected to have given continued elastic behavior. Thus we 
have /S.(j drop = o'i + — (Ti+i, where fi is the shear mod- 
ulus (determined as half the slope of the stress-strain curve), 
and AEdrop — Ei + ViaiAe — Ei^i, where Vi is the sys- 
tem volume. AEdrop (Aadrop) is positive if there is a drop 
in stress (energy). For some very small stress drops, the ap- 
parent energy drop is negative, due to finite resolution implied 
by a finite Ae. We count events with Addrop greater than a 
cutoff A(j,nin and AEdrop > 0. Our scaling analysis is based 
on the distributions of AEdrop and a quantity proportional to 
Aadrop^ that we call the slip volume, VgUp = V Aadrop /"^l-i- = 
V{Ae - (cr,+i - a,)/2fi) = V{Ae - Ae^;) = VAepi, where 
Aeei (Acpi) is elastic (plastic) strain. The significance of VsUp 
can be understood by supposing first that the plastic slip asso- 
ciated with an event is confined to a localized region of space, 
whose size had a narrow distribution independent of L. If the 
slip is characterized as a displacement d over an area A, then 
epi — Ad/V and VsUp = Ad. If an event involved m such 
elementary shear transformations the resulting VsUp would be 
^ niAd, and thus a measure of the number of elementary 
transformations which contribute to the macroscopic stress re- 
laxation. The cutoff Aajnin is chosen so that the minimum 
VsUp is independent of L and equal to SA'^. 

In Figs. |2] and [3] are shown cumulative distributions C{x) 
of AE and VsUp, where C{x) = P{x')dx', and P{x) is 
the probability distribution. The advantage of using C{x) is 
that it yields a much smoother curve, while no information 
is lost through binning. Furthermore, for power-law or ex- 
ponential behavior of P{x) at large x, C{x) maintains this 
behavior (with the exponent changing by one in the case of a 
power-law). The insets show P{x), which involves a binning 
procedure and result in noisy curves. The C{x) curves are al- 
most linear (with perhaps some downward curvature) in these 
semi-log plots, suggesting that the distributions are roughly 
exponential. The inverse-slopes, and hence the mean values, 
however, systematically increase with L, so they cannot be as- 
sociated with a characteristic event size independent of L. If 
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FIG. 2: (Color online) Left, cumulative probability distribution of 
A-Edrop in plastic events for different system sizes. Right, cumula- 
tive probability distribution of /S.Edrop/ L^'"^ ■ Inset, probability dis- 
tribution of A E drop- 
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FIG. 3: (Color online) Left, cumulative probability distribution of 
AVstip in plastic events for different system sizes. Right, cumulative 
probability distribution of AVsiip/ L^'^' ■ Inset, probability distribu- 
tion of AVsiip. 

these quantities are scaled by L^ '^ and L^-^, respectively, the 
C(x) curves collapse quite well onto master curves, as shown 
in the right panels of Figs. |2] and [3] Changing the exponent 
by 0.1 produces a slightly worse collapse in both cases. The 
means of AE and V^up scale as Li-43±o.03 ^^^^j ^i.63±0.04^ 
respectively; these exponents are consistent with the scaling 
collapses. 

Ref. jlH] interpreted the linear scaling with L in terms of 
the geometrical structure of the events: they tended to be ex- 
tended in one dimension, in the form of slip lines passing 
through the simulation cell. Extrapolating their results to 3D, 
one might expect planar events, scaling as L^. This is ex- 
cluded by our results. If the connection between the observed 
scaling of event distributions and the geometry of events is to 
be trusted, we can tentatively interpret the L"^/^ scaling as re- 
flecting a fractal geometry of the avalanches, somewhere be- 
tween string-like and planar Of course, the scaling of the 



form _L^/^ = y^/^ is very reminiscent of a central limit 
theorem — the variance in the extensive quantities E and Va 
should go like V (or N), and AE and Act correspond to the 
square root, the standard deviation. Note that the results of 
Ref. ll2lll are also consistent with this interpretation, since in 
two dimensions L is the square root of the system size. A 
central limit interpretation of this scaling would suggest that 
contributions from different parts of the system somehow add 
up in an uncorrected way. The clear change in event size dis- 
tribution with system size indicates that the events are spa- 
tially delocalized and in this way the resulting energy and 
stress drops certainly result from contributions from differ- 
ent parts of space. But within such a "random-noise" inter- 
pretation with uncorrected contributions from different parts 
of space one would expect that the atoms participating in the 
event would be more or less equally distributed over space. 
We shall see now that this is not the case. 

To study the geometrical structure of events, we need to 
define a measure of which atoms take part in an event. We 
choose to consider the atom displacements that take place dur- 
ing minimization rfj = |Afi|. Rather than impose an arbitrary 
cut-off to identify participating atoms, we use the participa- 
tion ratio P = {J2i df Y / {N ^ ■ d^), where N is the number 
of atoms. P = 1 for an event where all the di are equal and 
1/iV for one where a single atom moves; thus it is the effective 
fraction of participating atoms. In the left panel of Fig. |4]we 
show distributions of P for different L, as well as the means. 
The mean value of P is well fit by a power-law with exponent 
-1.44±0.03. To compare with the scaling of the extensive 
quantities AE and VgUp, we should consider VP = L^P, 
or add three to the exponent. This gives 1.56, close to the 
Vsiip exponent, implying a similar scaling for two measures 
of events which may both be considered "geometrical" in a 
sense — we saw above that VsUp is related to the amount of 
plastic strain at the boundaries that an event causes. 

Given P for an event, we define the set of set of partici- 
pating atoms as those whose di is in the top P of the popu- 
lation. This defines the "mobile" atoms without an arbitrary 
cut-off. We then define the fractal dimension Dp of the set in 
the usual box-counting way: Taking a box which contains all 
of the atoms, we divide the box by increasing powers of two 
(in all directions) and count how many boxes Nt for a given 
divisor d are required to contain all atoms in the subset. D p 
is determined as the best-fit slope on a plot of log{N\j) versus 
log{d), assuming the data lie on a straight line. Typically four 
points are available for the fit. The data generally exhibit a no- 
ticeable downward curvature, so the interpretation of a fractal 
geometry should not be taken too literally. Even so, it is in- 
teresting that the average D p computed this way indeed gives 
the value 1.6±0.05, independent of system size and close to 
the exponents determined from VsHp and P; distributions of 
Dp are shown in the right hand panel of Fig. H) These obser- 
vations support the idea that non-local plastic events can still 
be viewed as sums of similar elementary events, probably like 
elementary shear transformations. These sub-events organize 
in space along a fractal-like structure, and therefore must be 
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FIG. 4: (Color online) Left, distribution of participation ratios P for 
different system sizes L. Inset, the mean P versus L together with 
a power law fit with exponent -1.44. Right, distribution of apparent 
fractal dimension of participating atoms. 

strongly correlated. This suggests that mean-field approaches, 
such as the STZ theory[7], may be incomplete. 

In summary we have observed a clear signature of a ^3/2 
scaling of event sizes with system size, and presented evi- 
dence that this can be attributed to a fractal-like shape of the 
events. The range of system sizes is not large, but the statistics 
are good. It is interesting to speculate whether there is a rela- 
tion between the events seen in low temperature deformation 
and those in a supercooled liquid near the glass transition. The 
latter are known to exhibit strong spatial correlations; a recent 
attempt to infer a fractal dimension for individual clusters |j29[] 
yielded a value of 1.8 — not close enough to 1.5 to suggest an 
obvious connection, but enough perhaps to speculate that the 
dynamics changes in a smooth way upon going from high-T, 
zero stress to zero-T, high stress. Finally we mention the in- 
teresting question of how finite strain rates and temperatures 
cut off the scaling, and whether a length scale emerges from 
this cutting off, which could be connected to, for example, the 
observed width of shear bands (10-20 nm). 
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